library(plotly)
library(Hmisc)
library(WVPlots)
library(plyr)
source("tools/DataSplitTools.R")
source("tools/GeneralizedNadarayaWatson.R")
source("tools/CommonTools.R")
source("tools/CrossValTools.R")
Зчитаємо дані:
df <- read.csv2(file = "data/ZNOandVoating/input_2016.csv", header = 1, sep = ",", dec = ".", fileEncoding = "CP1251", stringsAsFactors = FALSE)
names(df) <- c("regname", "ukr","math", "pro_ukr", "radical", "opposition", "small", "not_voted")
head(df, 3)
## regname ukr math pro_ukr radical opposition small
## 1 Вінницька область 167 128 43.00824 6.005472 1.312608 7.75368
## 2 Вінницька область 137 0 43.00824 6.005472 1.312608 7.75368
## 3 Вінницька область 141 0 43.00824 6.005472 1.312608 7.75368
## not_voted
## 1 41.92
## 2 41.92
## 3 41.92
Одразу відмітимо, що в даних 1. Нема відомостей про голосування для Донецької та Луганської областей, а отже, вимушено не враховуємо абітурієнтів, що проживають у даних регіонах (майже 6% даних опуска’ється). 2. Київ викоремлено від Київської області
n_for_delete <- nrow(df[df$regname == "Донецька область" | df$regname == "Луганська область", ])
print(as.character(round(n_for_delete * 100 / nrow(df), 2) %&% "% записів, що відповідають Донецькій та Луганській облатям"))
## [1] "5.87% записів, що відповідають Донецькій та Луганській облатям"
df <- df[df$regname != "Донецька область" & df$regname != "Луганська область", ]
Виведемо загальні характеристики даних ЗНО з математики та української мови:
describe(df[, c("math", "ukr")])
## df[, c("math", "ukr")]
##
## 2 Variables 252262 Observations
## ---------------------------------------------------------------------------
## math
## n missing distinct Info Mean Gmd .05 .10
## 252262 0 55 0.77 53.53 69.96 0 0
## .25 .50 .75 .90 .95
## 0 0 123 157 172
##
## lowest : 0 100 104 107 111, highest: 196 197 198 199 200
## ---------------------------------------------------------------------------
## ukr
## n missing distinct Info Mean Gmd .05 .10
## 252262 0 83 0.997 125.4 58.72 0 0
## .25 .50 .75 .90 .95
## 111 137 165 184 189
##
## lowest : 0.0 100.0 101.0 103.0 104.0, highest: 197.0 198.0 199.0 199.5 200.0
## ---------------------------------------------------------------------------
Із таблиці вище видно, що в записах нема пропущених значень. Судячи із значень квантилів, половина значень ЗНО з математики є нулями. Варто відмітити, що значно менша частка значень 0 у ЗНО з української мови (близько 10%).
Подивимося на розподіли:
p1 <- plot_ly(x = df[, "math"], type = "histogram", name = "math")
p2 <- plot_ly(x = df[, "ukr"], type = "histogram", name = "ukr")
plotly::subplot(p1, p2) %>% layout(title = "EIT: Ukranian language and mathematics")
remove(p1); remove(p2)
Розподіли результатів ЗНО з української та математики значно різняться.
Подивимося на співвідношення людей, котрі отримали 0 по математиці та українській одночасно
print("Не подолали поріг або не здавали ЗНО відносно всіх абітурієнтів: " %&% as.character(round(nrow(df[(df$math == 0 & df$ukr == 0), ]) * 100 / nrow(df), 2)) %&% "%")
## [1] "Не подолали поріг або не здавали ЗНО відносно всіх абітурієнтів: 12.72%"
print("Не подолали поріг або не здавали ЗНО відносно тих, хто не склав математику: " %&% as.character(round(nrow(df[(df$math == 0 & df$ukr == 0), ]) * 100 / nrow(df[(df$math == 0), ]), 2)) %&% "%")
## [1] "Не подолали поріг або не здавали ЗНО відносно тих, хто не склав математику: 20.77%"
print("Не подолали поріг або не здавали ЗНО відносно тих, хто не склав українську: " %&% as.character(round(nrow(df[(df$math == 0 & df$ukr == 0), ]) * 100 / nrow(df[(df$ukr == 0), ]), 2)) %&% "%")
## [1] "Не подолали поріг або не здавали ЗНО відносно тих, хто не склав українську: 91.19%"
Із останнього напрошуєтся висновок, що ті, хто не склали іспит із української, скоріш за все не склали й математику. Що є надзвичайно логічно: для вступу в університет абітурієнту необхідно отримати певний мінімум балів із кожного ЗНО. Першою завжди здається українська. Тому, не склавши українську, сенс складати математику залишається у винятковим.
Поглянемо на кореляції між результатами ЗНО, політичними силами:
corrolation_m <- cor(df[, -1], method = "pearson")
plot_ly(x = rownames(corrolation_m), y = colnames(corrolation_m), z = corrolation_m, colors = "Greys", type = "heatmap") %>%
layout(title = "Pearson Correlation")
# remove(corrolation_m)
Відмітимо майже 0 кореляцію між результатами ЗНО з української мови та про-українським політичним напрямком та слабку кореляцію між результатми математики та української.
Досить висока кореляція між результатми не голосування та голосуванням за малі партії. Така ситуація може пояснюватися тим, що населення не підтримує основних лідерів, а тому не йде голосувати або голосує спеціально за тих, хто не пройде.
Видалимо стрічки, що відповідають спостереженням із 0 хочаб в одному ЗНО та подивимося на розподіл.
df <- df[(df$math != 0) & (df$ukr != 0) ,]
p1 <- plot_ly(x = df[, "math"], type = "histogram", name = "math")
p2 <- plot_ly(x = df[, "ukr"], type = "histogram", name = "ukr")
plotly::subplot(p1, p2) %>% layout(title = "Distribution of the passed EIT",
annotations = list(list(x = 0.2, y = 1.05, text = "math", showarrow = F, xref = 'paper', yref = 'paper'),
list(x = 0.8, y = 1.05, text = "ukr", showarrow = F, xref = 'paper', yref = 'paper')))
remove(p1); remove(p2);
Відмітимо, що наразі дуже відсіялася частка людей, котрі погано здали українську мову (було в околі 3-х тисяч, стало в около 700).
Гістограма розсіювання, в свою чергу, виглядає наступним чином:
plot_ly(data = df, x = ~math, y = ~ukr, mode = "markers", type = "scatter", marker=list(opacity = 0.7)) %>% layout(title = "Passed EIT results")
Варто відмітити, що для абітурієнтів, що склали математику краще ніж 160 балів, помітна тенденція, як правило, кращого складання ЗНО з математики.
Подивимося на розподіли по кожній політичній силі окремо (якщо конкретно обрана політична сила має для даного спостереження найбільший вплив (найбільшу ймовірність), то дане спостереження відносимо до ціє політичної сили)
W <- as.data.frame(df[, -c(1, 2, 3)])
W[, "max_value"] <- apply(W, 1, max)
all_plots <- list()
annotations_list <- list()
names_ <- names(W[, -ncol(W)])
for (i in seq_along(names_)) {
name <- names_[i]
d <- df[W[, name] == W["max_value"], ]
if (nrow(d) != 0) {
all_plots[[name]] <- plot_ly(x = d[, "math"], type = "histogram", name = name)
annotations_list[[i]] <- list(x = 0.15*i, y = 1.005, text = name, showarrow = F, xref = 'paper', yref = 'paper')
}
}
plotly::subplot(all_plots) %>%
layout(title = "Mathematics EIT according to political strength", annotations = annotations_list)
remove(d); remove(all_plots); remove(annotations_list); remove(i); remove(name)
Відмітимо, що візуально розподіли балів виглядають майже однаково для про-украхнських регіонів та тих, де в основному, не голосували.
all_plots <- list()
annotations_list <- list()
for (i in seq_along(names_)) {
name <- names_[i]
d <- df[W[, name] == W["max_value"], ]
if (nrow(d) != 0) {
all_plots[[name]] <- plot_ly(x = d[, "ukr"], type = "histogram", name = name)
annotations_list[[i]] <- list(x = 0.15*i, y = 1.005, text = name, showarrow = F, xref = 'paper', yref = 'paper')
}
}
plotly::subplot(all_plots[["pro_ukr"]], all_plots[["not_voted"]]) %>%
layout(title = "Ukr EIT according to political strength", annotations = annotations_list)
remove(d); remove(all_plots); remove(annotations_list); remove(names_); remove(name); remove(i)
Розподіли результатів ЗНО з уркаїнської візуально більш різняця для різних політичних вподобань.
for (name in names(W[, -ncol(W)])) {
d <- df[W[, name] == W["max_value"], ]
if (nrow(d) != 0) {WVPlots::ScatterHist(d, "math", "ukr", name, smoothmethod = "none", estimate_sig = FALSE, minimal_labels = TRUE)}
}
remove(d); remove(W); remove(name); remove(df)
Із останньої пари зображень здається, що нижній поріг результатів української та математики (для абітурієнтів із математикою > 160 балів) вищий у випадку про-українських регіонів. Така ситуація можлива також через різну кількість абітурієнтів у порівнювальних регіонах.
Виведемо загальні характеристики кожної політичної сили:
df <- read.csv2(file = "data/ZNOandVoating/input_2016.csv", header = 1, sep = ",", dec = ".", fileEncoding = "CP1251", stringsAsFactors = FALSE)
names(df) <- c("regname", "ukr","math", "pro_ukr", "radical", "opposition", "small", "not_voted")
df <- df[df$regname != "Донецька область" & df$regname != "Луганська область", ]
describe(df[, -(1:3)])
## df[, -(1:3)]
##
## 5 Variables 252262 Observations
## ---------------------------------------------------------------------------
## pro_ukr
## n missing distinct Info Mean Gmd .05 .10
## 252262 0 23 0.997 33.51 13.31 15.71 16.23
## .25 .50 .75 .90 .95
## 21.52 34.08 42.46 48.82 53.87
##
## lowest : 15.71244 16.22691 19.37149 20.13312 21.00261
## highest: 43.00824 44.77893 48.82355 50.41795 53.87200
## ---------------------------------------------------------------------------
## radical
## n missing distinct Info Mean Gmd .05 .10
## 252262 0 23 0.997 6.703 2.561 3.047 3.832
## .25 .50 .75 .90 .95
## 4.992 6.580 8.477 8.912 10.428
##
## lowest : 3.046992 3.831828 3.834072 4.523960 4.619912
## highest: 8.886996 8.911950 10.084956 10.427880 11.555155
## ---------------------------------------------------------------------------
## opposition
## n missing distinct Info Mean Gmd .05 .10
## 252262 0 23 0.997 4.295 4.541 0.4165 0.4970
## .25 .50 .75 .90 .95
## 1.3126 2.0612 7.1334 11.6156 14.5749
##
## lowest : 0.344142 0.416508 0.497000 0.953295 1.042825
## highest: 6.796640 7.133360 10.466742 11.615622 14.574912
## ---------------------------------------------------------------------------
## small
## n missing distinct Info Mean Gmd .05 .10
## 252262 0 23 0.997 9.319 2.148 6.779 7.189
## .25 .50 .75 .90 .95
## 7.754 9.167 9.966 13.113 13.113
##
## lowest : 5.856787 6.778902 7.189344 7.360584 7.496422
## highest: 10.890088 11.198572 11.346280 13.112736 13.519935
## ---------------------------------------------------------------------------
## not_voted
## n missing distinct Info Mean Gmd .05 .10
## 252262 0 23 0.997 46.17 9.604 30.00 35.15
## .25 .50 .75 .90 .95
## 41.92 44.14 52.81 57.20 60.48
##
## lowest : 30.00 31.72 35.15 36.27 39.79, highest: 54.68 55.32 57.20 58.64 60.48
## ---------------------------------------------------------------------------
Одразу видно, що 23 унікальні значення в кожному стовпчику відповідають 22 областям України та м. Києву. Варто відмітити, що відсутні області, із переважними радикальними або опозиційними силами. Також немає регіонів, у яких голоса віддані виключно за малі політичні партії.
Розподіл абітурієнтів по областях:
ab <- plyr::ddply(.data = df, .variables = .(pro_ukr, radical, opposition, small, not_voted), .fun = nrow)
plot_ly(ab, y = ~V1, type = "bar") %>%
layout(title = "Abiturients distribution according regions",
xaxis = list(title = "Region ID"), yaxis = list(title = ""))
Варто відмітити, що кілька регіонів надто відрізняються чисельністю. Через таку зміщенність даних можуть бути зміщені результати.
Із розподілу відсотків підтримки кожної політичної сили відповідно до регіону видно, що більшість регіонів є політично неактивними. Також варто відмітити, що розподіли голосів малих партії, у більшості випадків, збігається з розподілом голосів радикальної партії.
plot_ly(y = ab[, "pro_ukr"], type = "bar", name = "pro_ukr") %>%
add_trace(y = ab[, "not_voted"], name = "not_voted") %>%
add_trace(y = ab[, "radical"], name = "radical") %>%
add_trace(y = ab[, "opposition"], name = "opposition") %>%
add_trace(y = ab[, "small"], name = "small")
remove(ab); remove(df)
На основі порівняння розподілів можна зробити висновок про те, що дані можуть бути описані за допомогою моделі суміші зі змінними концентраціями. Проте, зважаючи на кореляції та відсоткові співвідношення, можливо, модель суміші зі змінними концентраціями, буде не надто ефективною.
Розіб’ємо вибірку на тренувальну та тестову частини у відсотковому співвідношенні 80/20%. Далі, тренувальну частину розіб’ємо на 5 рівних частин для проведення 5-ти фолдової кроссвалідації з метою визначення оптимального для кожної компоненти параметра згладжування h. Узагалі кажучи, можливо, варто було б врахувати незбалансованість даних під час розбиття, проте ми цього не робили.
Під час 2-х компіляцій даного звіту для різних розбиттів на тренувальні-тестові дані, було отримано різні якості прогнозування побудованих моделей.
Для того, щоб порівняти як у цілому справляється запропонований метод прогнозування на даних, що наразі досліджуються, проведемо наступний експеримент:
Отже, розіб’ємо на тренувальну-тестові вибірки та збережемо їх.
set.seed(42)
df <- read.csv2(file = "data/ZNOandVoating/input_2016.csv", header = 1, sep = ",", dec = ".", fileEncoding = "CP1251", stringsAsFactors = FALSE)
names(df) <- c("regname", "ukr","math", "pro_ukr", "radical", "opposition", "small", "not_voted")
df <- df[df$regname != "Донецька область" & df$regname != "Луганська область", ]
df <- df[, -1]
df[, -(1:2)] <- df[, -(1:2)]/100
df <- df[(df$math != 0) & (df$ukr != 0) ,]
train_test_list <- train_test_split(df, ratio = 0.80)
train <- train_test_list[["train"]]; test <- train_test_list[["test"]]
write.csv(train, "data/computation_results/train.csv"); write.csv(test, "data/computation_results/test.csv")
remove(df); remove(train_test_list); remove(test)
Скористаємося 5-ти фолдовою кросс-валідацією для знаходження оптимальних параметрів згладжування для кожної моделі (проведемо описаний нижче есперимент 5 разів).
Для визначення оптимального параметра використаємо модуль МНК зважений ваговими коефіцієнтами:
\[MSE^{(m)} = \frac{1}{N} | \sum_{i=1}^{N}{a_{i:N}^{(m)}(Y_i - \hat{g}^{(m)}(X_i))^2} |\] Зауважимо, що модуль використовується для того щоб підчас підрахунків середнього зваженого MSE по всіх фолдах фіксованої компоненти, фіксованого h, не отримати 0 (який може бути хибно протрактований).
Отже, визначення оптимальних параметрів h відбувається за схемою описаною нижче.
Зафіксуємо h, порахуємо зважений МНК для кожного фолда та знайдемо середнє та дисперсію по 5-ти фолдам для кожної компоненти. (Так прорахуємо для кожного фіксованого h.)
Відкинемо 25% значень середніх, які відповідають найбльшим значенням дисперсії. Таким чином залишаємо лише ті результати, на яких моделі показали більш стабільну поведінку.
Серед значень середніх результатів моделей визначаємо h, які відповідають найменшим середнім зваженим МНК (окремо для кожної компоненти).
# h_range <- c(seq(0.1, 1, 0.2), seq(1, 5, 0.5), 8, 10)
#
# for (i in 1:5) {
# print("***** " %&% as.character(i) %&% " ******")
# set.seed(i)
# cv_split <- cross_validation_split(train, k = 5)
# GNW_cv_results <- GNWcv_across_h(h_range = h_range,
# cv_df_split = cv_split,
# X_colname = "math", Y_colname = "ukr",
# W_colname = c("pro_ukr", "radical", "opposition", "small", "not_voted"),
# use_parallel = TRUE)
# write.csv(GNW_cv_results, "data/computation_results/2016cc_data/5cv_mean_std_" %&% as.character(i) %&% ".csv")
# }
#
# # h_opt <- optimal_h(GNW_cv_results)
# remove(h_range); remove(i); remove(cv_split); remove(GNW_cv_results)
Результат середніх та диспресії збережемно в змінній \(GNW\_cv\_results\), а оптимальні \(h\), серед запропонованих для перебору, у змінну \(h\_opt\).
data_path <- "data/computation_results/2016cc_data/"
file_names <- dir(data_path) #where you have your files
file_names <- data_path %&% file_names
cv_results <- lapply(file_names,
function(x){
df <- as.matrix(read.csv(x, row.names = 1))
l <- list(df=df , h=optimal_h(df))
return(l)
})
cv_results_h <- lapply(cv_results, function(x){x$h})
cv_results_mean_std <- lapply(cv_results, function(x){x$df})
print(cv_results_h)
## [[1]]
## [,1] [,2] [,3] [,4] [,5]
## [1,] 438.4071 454.3044 1636.697 1035.322 351.3797
## [2,] 2.5000 10.0000 0.100 5.000 2.5000
##
## [[2]]
## [,1] [,2] [,3] [,4] [,5]
## [1,] 752.0023 610.5381 1900.39 522.0256 1014.304
## [2,] 5.0000 5.0000 10.00 8.0000 4.000
##
## [[3]]
## [,1] [,2] [,3] [,4] [,5]
## [1,] 636.954 574.1809 770.329 409.463 539.1719
## [2,] 4.500 10.0000 0.100 3.500 0.1000
##
## [[4]]
## [,1] [,2] [,3] [,4] [,5]
## [1,] 261.346 488.8267 1669.38 476.3409 145.2318
## [2,] 3.000 5.0000 2.50 4.5000 4.0000
##
## [[5]]
## [,1] [,2] [,3] [,4] [,5]
## [1,] 485.6412 512.2064 1866.826 271.5143 205.7276
## [2,] 3.5000 5.0000 5.000 8.0000 2.5000
Порівняємо тепер результати кросс-валідації з результатми прогнозування на тестовій частині даних.
# train <- read.csv("data/computation_results/train.csv", row.names = 1, stringsAsFactors = FALSE)
# test <- read.csv("data/computation_results/test.csv", row.names = 1, stringsAsFactors = FALSE)
# h_opt <- cv_results_h[[4]][2, ]
#
# GNW <- GeneralisedNadarayaWatson$new()
# GNW$train(X_train = as.numeric(train[, "math"]),
# Y_train = as.numeric(train[, "ukr"]),
# W_train = as.matrix(train[, -c(1, 2)]))
# GNW$run_cluster()
# prediction_with_acoeff <- GNW$predict_in_parallel(X_test = as.numeric(test[, "math"]),
# W_test = as.matrix(test[, -c(1, 2)]), h = h_opt)
# GNW$stop_cluster()
#
# prediction <- prediction_with_acoeff[["prediction"]]
# A_coeff <- prediction_with_acoeff[["A_test"]]
#
# write.csv(prediction, "data/computation_results/prediction2016v.csv")
# write.csv(A_coeff, "data/computation_results/A_coeff2016v.csv")
#
# remove(GNW); remove(train)
Порахуємо зважену суму МНК для прогнозу зробленого на тестових даних.
train <- read.csv("data/computation_results/train.csv", row.names = 1, stringsAsFactors = FALSE)
test <- read.csv("data/computation_results/test.csv", row.names = 1, stringsAsFactors = FALSE)
prediction <- as.matrix(read.csv("data/computation_results/prediction2016v.csv", row.names = 1, stringsAsFactors = FALSE))
A_coeff <- as.matrix(read.csv("data/computation_results/A_coeff2016v.csv", row.names = 1, stringsAsFactors = FALSE))
WMSE <- weighted_MSE(Y_true = as.numeric(test[, "ukr"]), Y_predicted = prediction, A_coeff = A_coeff)
names(WMSE) <- rep("comp.") %&% as.character(1:5)
WMSE
## comp.1 comp.2 comp.3 comp.4 comp.5
## 517.8019 -449.8565 366616.4095 -553.3369 -1092.8865
Відмітимо, що найменше значення зваженого МНК для моделей, що відповідають проукраїнській силі. Нагадаємо, найбільша по модулю кореляція спостерігалася між результатами ЗНО з української мовои та проукр. силами й тими, хто не проголосував. Також відмітимо, що раніше ми змогли виділити лише регіони, де переважає проукр напрямок та ті, в яких переважає не голосування. Тобто не виявилось регіонів, у яких переважають голоси за радикалів, за опозицію або за малі партії.
Порівняємо зважені МНК результатів прогнозування з результатами крос-валідації, використовуючи відповідними параметрами згладжування.
h_opt <- cv_results_h[[4]][2, ]
GNW_cv_results <-cv_results_mean_std[[4]]
cv_mean_results <- sapply(1:length(h_opt),
function(i){
GNW_cv_results[which("mean_" %&% as.character(h_opt[i]) == rownames(GNW_cv_results)), i]
})
cv_test_comparing <- rbind(cv_mean_results, WMSE)
rownames(cv_test_comparing) <- c("cross-validation", "test")
colnames(cv_test_comparing) <- 1:5
print(cv_test_comparing)
## 1 2 3 4 5
## cross-validation 261.3460 488.8267 1669.38 476.3409 145.2318
## test 517.8019 -449.8565 366616.41 -553.3369 -1092.8865
Із останньої таблиці видно, що квадрати залишків першої моделі, порахованих для кросс-валідації та на відкладеній тестовій вибірці не надто різняця у порівнянні з іншими компонентами.
Це свідчить про більш-менш вдалий підбір параметрів згладжування для моделі, що відповідає проукраїнському настрою.
Для 3 інших моделей результати - невтішні. Останнє може бути зумовлене слабкою представленістю відповідних класів у вибірці. Можливо є сенс розглянути 3-х компонентну суміш: проукраїнську направленість, “не проукраїнську” (об’єднати результати голосуваня за малі партії, опозицію, радикалів) та “байдужих”.
Зобразмо спочатку загальну картину: гістограму розсіювання з усіма прознозами для всіх моделей:
plot_ly(data = test, x = ~math, y = ~ukr, mode = "markers", type = "scatter", name="EIT results") %>%
add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 1], name = "pro_ukr", mode = 'lines+markers') %>%
add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 2], name = "radical", mode = 'lines+markers') %>%
add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 3], name = "opposition", mode = 'lines+markers') %>%
add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 4], name = "small", line=list(color="yellow"), marker=list(color="yellow", line=list(color="yellow")), mode = 'lines+markers') %>%
add_trace(x=sort(test[, "math"]), y = prediction[sort(test[, "math"], index.return = TRUE)$ix, 5], name = "not_voted", mode = 'lines+markers') %>%
layout(title = "Passed EIT results (test dataset)")
Окремо зобразимо для абітурієнтів, які проживають в областях, що ймовірніше наслідують проукраїнську групу:
train <- read.csv("data/computation_results/train.csv", row.names = 1, stringsAsFactors = FALSE)
test <- read.csv("data/computation_results/test.csv", row.names = 1, stringsAsFactors = FALSE)
prediction <- as.matrix(read.csv("data/computation_results/prediction2016v.csv", row.names = 1, stringsAsFactors = FALSE))
A_coeff <- as.matrix(read.csv("data/computation_results/A_coeff2016v.csv", row.names = 1, stringsAsFactors = FALSE))
Wtest <- as.data.frame(cbind(prediction, test[, -c(1, 2)]))
Wtest[, "max_value"] <- apply(Wtest[, -(1:5)], 1, max)
pro_ukr <- test[Wtest[, "pro_ukr"] == Wtest["max_value"], ]
p1 <- plot_ly(data = pro_ukr, x = ~math, y = ~ukr, mode = "markers", type = "scatter", name = "EIT results in pro_ukr regions") %>%
add_trace(x = test[rownames(pro_ukr), "math"],
y = Wtest[rownames(pro_ukr), 1], name = "pro_ukr model", mode = 'markers')
p1
Відмітимо, що основні тренди моделі для проукраїнських регіонів не надто збереглися.